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Abstract. We study the fluctuations of the current J(t) of the totally asymmetric 
exclusion process with open boundaries. Using a density matrix renormalization group 
approach, we calculate the cumulant generating function of the current. This function 
can be interpreted as a free energy for an ensemble in which histories are weighted by 
exp(— sJ(t)). We show that in this ensemble the model has a first order space-time 
phase transition at s = 0. We numerically determine the finite size scaling of the 
cumulant generating function near this phase transition, both in the non-equilibrium 
steady state and for large times. 
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1. Introduction 

Heat or particle currents arise in systems that are driven away from equilibrium by 
bringing them in contact with reservoirs at different temperatures or chemical potentials. 
In macroscopic systems, fluctuations of these currents can often be neglected and a 
description using non-equilibrium thermodynamics is then appropriate [1J. However, 
current fluctuations can become important in mesoscopic systems and in the vicinity of 
non-equilibrium critical points. 

In recent years, the interest in current fluctuations has increased considerably. 
From the experimental side, it is nowadays possible to measure third- and higher- 
order cumulants of the current in problems of charge transport [21 E]. In quantum 
point contacts, current fluctuations can be used as an entanglement meter [4]. 
Theoretically, one is interested in the entropy produced in the presence of a current, 
since its fluctuations have been shown to obey various kinds of fluctuation theorems 

ei m ei ei m Ell- 
in this article, we study current fluctuations in the totally asymmetric exclusion 
process (TASEP). In this classical, stochastic model, particles on a one-dimensional 
lattice jump to the right (provided that site is empty), thus giving rise to a time- 
dependent particle current. The TASEP was originally introduced to describe the 
biological process of translation, in which ribosomes move along a messenger RNA 
and use the information stored on it to build proteins [12]. In this context, one can 
relate current fluctuations to fluctuations in protein production as we recently showed 
[13j . Nowadays, the TASEP, and some of its variants, have become paradigmatic for 
non-equilibrium statistical mechanics [14]. 

Fluctuations in the current can be determined from the cumulant generating 
function (CGF). In the large time limit, the TASEP evolves to a non-equilibrium steady 
state (NESS) for which the CGF equals the largest eigenvalue of a generalized generator 
[14] . This largest eigenvalue has been determined exactly for the TASEP on a ring 
[T5] [To] using the Bethe-ansatz and was found to have an interesting scaling form, in 
which the KPZ-dynamical exponent z = 3/2 [T71 IT8] appears. Much less is known for 
a configuration with open boundaries where particles can enter on the left, and leave 
on the right. The average and variance of the current have been determined using a 
matrix approach [191 EQ], but no results are known for higher cumulants or for the full 
generating function. Recently, we showed how the latter can be determined numerically 
using a density matrix renormalization group (DMRG) approach |2T] . In this article, 
we further analyse our results and show that in the thermodynamic limit, the CGF has 
a non-analytic behaviour which corresponds to a space-time phase transition of a type 
that has been discovered recently in models of glassy dynamics [22l |23| [24] [25] . We 
establish the scaling form of the CGF near this phase transition. 

Using the same DMRG approach we also investigate the gap in the spectrum of 
the generalised generator. This gap determines the approach of the current fluctuations 
to their NESS-value. We investigate the scaling properties of the gap and conjecture a 
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scaling form for the time dependent current fluctuations which so far have been hardly 
studied, but which could be relevant, for example, in applications to protein production 
|13j . The scaling form is verified with simulations based on the Gillespie algorithm. 

This paper is organized as follows. In section 2 we introduce the model, the 
cumulant generating function of the current, and the reinterpretion of this function 
as a free energy in the so called s-ensemble. In section 3, we show that the CGF 
must be a non-analytical function of the variable s. In section 4 we discuss the DMRG 
approach. In the section 5 and 6 we present the results of our numerical calculations 
on the CGF and the time dependence of the current fluctuations respectively. We pay 
particular attention to the scaling behaviour of these quantities. Finally, in section 7 
we present our conclusions. 



2. The cumulant generating function and the s-ensemble for the TASEP 

In the totally asymmetric exclusion process, each site of a one-dimensional lattice of L 
sites can be empty or occupied by at most one particle. The dynamics of the TASEP 
is that of a continuous time Markov process for which the probability P(C, t) that the 
system is in a configuration C at time t evolves according to the master equation 

d t P(C, t) = E [w(C -+ C)P(C, t) - w(C -> C')P(C, t)} (1) 

Here w(C — > C) is the transition rate from configuration C to C In the TASEP, 
particles can only jump to the right (provided that site is empty) with unit rate. 
We impose open boundary conditions for which particles can enter the lattice on the 
left-side with rate a < 1, and leave it on the right with rate < 1. For further 
reference, we introduce the inverse lifetime, or escape rate, r(C) of the state C, 
r(C) = Ec^cw(C^C). 

As a consequence of these dynamics, a current flows from left to right through the 
system. The current at time t, j{i,t), gives the number of particles passing per unit of 
time through the bond between sites % and i + In a particular realisation, or history, 
of the TASEP, we can count the total integrated current J(t, L) = ^2 i=0 Joj(i, t')dt' > 
through all bonds up to time t (where the particles entering respectively leaving the 
system correspond with i — and i = L). This is a stochastic variable whose properties 
can be obtained from the cumulant generating function 

v(s,t,L) = hn(e- sJ ^) (2) 

where the average (•) is taken over the histories of the process. From ([2]), the average 
current per unit time J(t, L) [fj its variance AJ(t, L) and higher cumulants can be found 
by taking derivatives: 

7{tj7) =i(J(t,L)) = -^(0,t,L) (3) 

\ overlined quantities represent averages over time, quantities with a * represent time averages in the 
NESS. 
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AJ(t,L) = I \(J\t,L)) - (J(t,L)f] = P^(0,t,L) (4) 



t 



ds 2 



Equation ^ has a large mathematical similarity with the definition of the dimensionless 
free energy in equilibrium situations in which the inverse temperature (3 = 1/ksT is 
replaced by s and the sum over microstates is replaced by a sum over histories. We 
can therefore interpret Q as a thermodynamics of histories [26] in which realisations 
of the stochastic process are weighted by exp (— sJ(t, L)). For negative s, ^ will be 
dominated by histories with a large current, whereas for positive s it will mainly get 
contributions from those with a small current. There is no clear physical meaning to the 
variable s, but in the context of glassy dynamics, it has been shown that extending the 
parameter space with the variable s can lead to new and interesting insights [22l 1231 [24] . 

In the long time limit, the TASEP evolves to a unique non-equilibrium steady 
state (NESS). The average current J*(L) = lim^oo J(t, L) in this state can be obtained 
exactly with a matrix technique [12]. If one also takes the thermodynamic limit, the 
current per bond J* = Hindoo J*(L)/(L + 1) is found to behave non-analytically as a 
function of a and (3 leading to the recognition of three phases [19]. In the low-density 
(LD) phase (a < 1/2,(3 > a), J* = a(l — a), while in the high-density (HD) phase 
((3 < 1/2, a > (3), J* equals (3(1 — (3). Finally, in the maximal current (MC) phase 
(a > 1/2,(3 > 1/2), J* = 1/4. The TASEP therefore has three (boundary driven) 
nonequilibrium phase transitions. 

It is useful to also introduce the s-weighted average current in the NESS together 
with its variance which are defined as 

r(s,L) = -\im^(s,t,L) (5) 

A*J(s,L)= \im ^(s,t,L) (6) 

For clarity, we now repeat a standard argument [26] that shows that in the NESS, 
the cumulant generating function equals the largest eigenvalue of a matrix H(s), and 
that the approach to the asymptotic value is determined by the lowest gap of that 
matrix. We therefore introduce firstly the probability P(C, J,t) that the system is in 
configuration C and has integrated current J at time t. Using 0, we immediately find 
that 

d t P(C,J,t)= ]T [w(C ^C)P(C',J -l,t)-w(C ^C')P(C,J,t)}(7) 

Consequently, the discrete Laplace transform P(C,s,t) = jLo e' sJ P(C, J, t) evolves 
according to 

d t P{C,s,t)= [ W ( C ' ->• C)e- s P{C',s,t) -w{C ->• C')P(C,s,t)] (8) 

C'+C 

To continue, it is convenient to introduce a matrix notation as is common in the so 
called 'quantum' approach to stochastic particle systems [27J. We therefore introduce a 
set of basis vectors \C) each corresponding to a configuration C and a vector \P(s,t)) 
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with components P(C,s,t) = (C\P(s,t)}. Using this notation, the set of equations (|8 
can be rewritten as 

d t \P(s,t)) =H(s)\P(s,t)) (9) 

The diagonal elements of the matrix H(s) are equal to minus the inverse lifetimes of 
the states, while the off-diagonal elements are given by the transition rates multiplied 
by e~ s . For s — 0, d9j) reduces to the master equation (jl]) and H(0) corresponds to the 
generator of the stochastic process. We will therefore refer to H(s) as the generalized 
generator. 

Using the 'quantum' notation of [27] in which an empty (occupied) site is associated 
with an up (down) spin, one can easily show that for the TASEP 

L-l 



+ oi 



e s s 1 



(10) 



i=i 



Here ni : v i: sf and Sj are standard particle number, vacancy number, particle 
annihilation and creation operators at site i 



n 




1 



1 




1 






1 



The formal solution to ^ is 

\P(s,t)) = e H ^\P(s,0)) 
Therefore, we have 

( e -J(*,i)) = J2Y, e ~ Sj(t,L)p ( C > 
c J 

c 

= Y,(C\e H{s)t \P(s,0)) 

c 



(12) 



(13) 



Using the spectral theorem, the sum in (13) can be written in terms of the eigenvalues 
and eigenvectors of H(s). In the long time limit, this sum is dominated by the 
largest eigenvalue Ai(s, L) while the first correction term involves the gap G(s,L) = 
Xi(s,L) — X 2 (s,L) with the second largest eigenvalue X 2 (s,L). One has 



( e -»J(*,i)) = A x e x ^ L)t \l + A 2 e 



G(s,L)t 



+ 



(14) 



where A\ and A 2 are time-independent factors depending on the initial conditions and 
the eigenvectors associated with the two largest eigenvalues. Therefore, for t — > oo 



At*(s,L) = ]imft(s,t,L) = lim - ln( e - sJ( *' L) ) = Ai(s,L) 



(15) 



which shows that fi* is intensive in time. 
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3. Space-time phase transition 

In this section, we give a simple argument that shows that in the NESS and in the 
thermodynamic limit the cumulant generating function is non-analytic at s — 0. Such 
phase transitions in the properties of histories have been called space-time transitions. 

Firstly, we observe that fi(s, t, L) is a non-increasing function of s and by definition 
is zero at s = 0. Hence we obtain the bounds 

> fi(s,t,L) > n(s ->■ oo,t,L) s>0 (16) 

For s oo, only histories with J(t, L) = contribute to /i. This implies that the system 
at time t is still in the configuration Co in which it was initially. The probability for 
this to happen decays exponentially with waiting time l/r(Co). Taking an average over 
possible initial conditions, one obtains 

Jim ( e -W)> = Y,Po(C )P(C ,J = 0,t) = E^oK^' (17) 

Co Co 

where po(Co) is the probability that the system is at t = in Cq. In the long time limit, 
this sum will be determined by the configurations with the largest lifetime. It is not too 
difficult to realise that these are the completely empty configuration C e with r(C e ) = a 
and the fully occupied configuration Cj with r(Cf) = (3. Therefore, one has for very 
large times 



lim (e- sJ ^)= Pmm e- mi ^ t 



I _j_ Pmax c -\ a -p\t + 
Pmin 



(where p min = p (C e ) and p max = p (C f ) if a < (3, and p max = p (C e ) and p min = p (C f ) 



if a > /3.) Comparing with (14) we have that Ai(s — >■ oo, L) = — min[a,/3] and 



G(s — > oo,L) = \a — f3\. Inserting (18) in (16) and taking t — > oo, one obtains 



> fi*(s, L) > - min[a, (3} (19) 

These inequalities imply that the space-intensive quantities lim^oo /i*(s, L)/(L + 1) 
and lim^oo J* (s,L)/(L + 1) are equal to zero for any strictly positive s. On the 
other hand, as discussed above, it is known from the exact solution [19] that J* = 
linit^oo J* (0, L)/(L + 1) is non-zero. Hence, it follows that in the s-ensemble, the 
TASEP has a first order space-time phase transition at s — 0, and this for every a and 
P. 

It is interesting to remark that this transition may be absent in the partially 
asymmetric exclusion process. In that model, particles can also jump to the left with 
a rate q and enter the system on the right and leave at the left side. Hence no lower 
bound on J(t, L) can be given and the argument presented above does not hold. 

In the remainder of this paper, we are interested in the finite size scaling properties 
of the cumulant generating function near the transition at s = 0. For the TASEP on 
a ring it has been shown exactly (T5J [16] that the cumulant generating function in the 
NESS scales as 



ti ing (s, L) = sLp(l - p) + W ^-^#(^27^(1 - p)L») (20) 



Current fluctuations in the TASEP 



7 



Here p is the density of particles and if is a scaling function. The first term on the 
right hand side of this equation equals sLJ* since in the ring case J* — p(l — p). In 
the remainder of this paper, we will use a numerical approach to investigate the scaling 
properties of p for the case of open boundaries. 



4. The DMRG approach 

The behaviour of the cumulant generating function of the current at large times in a 
finite system is determined by the two largest eigenvalues of the generalized generator 



(10). If we interpret this generator as a 'Hamiltonian' [27], we realise that calculating 
these eigenvalues is mathematically similar to determining the ground state energy and 
the gap of a quantum spin chain. The main difference from a standard quantum problem 
comes from the non-hermiticity of the generalised generator. 

One of the most precise numerical procedures for calculating ground state properties 
of quantum chains is the density matrix renormalization group (DMRG) introduced by 
S White [2HI ESI ED]. Later, this technique was adapted to generators of stochastic 
processes [31] 132] and in this context a study of the gap in the generator (at s — 0) of 
the TASEP was made [33J. More recently, we applied the DMRG for the first time to 
the generalized generators associated with current fluctuations in driven lattice gases 
and activity fluctuations in the contact process [2TJ US] • 

The application of the DMRG to these generalized generators is not fundamentally 
different from the standard approach used for quantum systems. For some of the 
technical aspects, we refer to [31]. 

We have used the DMRG to calculate the two largest eigenvalues of the generator 



(10) as a function of s for various points in the phase diagram of the TASEP. We are 
typically able to reach reliable results for L up to 60. At s = it is possible to go up to 
L « 100. The upper limit for L that can be reached is essentially set by the stability of 
the Arnoldi algorithm used to diagonalise large non-Hermitian matrices. 

Once the largest eigenvalue have been calculated with sufficient numerical accuracy, 
the average current and its variance in the s-ensemble are determined using numerical 
differentiation. Because of numerical round-offs, it is not possible to obtain cumulants 
beyond the second with sufficient accuracy. 



5. Current fluctuations in the stationary state 

In this section, we present our results for the current fluctuations in the NESS. 

In Figure [I] we show a typical result for the CGF p*(s,L) = Xi(s,L). The data 
are for a point in the low density phase with a = 0.35, (3 = 2/3. In the inset, it can be 
seen that already for rather small s-values the asymptotic value — min[a, 0] = —0.35 is 
reached, as predicted in the previous section. When L increases this limiting value is 
reached for ever smaller s. The behavior shown here is typical for all values of a and f3 
that we investigated. 
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Figure 1. (Color online) Cumulant generating function as a function of s for 
a = 0.35,(3 — 2/3 for L = 10,20,30,40,50 and 60 (increasing from bottom to top 
for s < ). The region with < s < 1 is shown in more detail in the inset (L — 20, 40 
and 60, increasing from top to bottom). On the right we show the CGF per bond for 
the same system sizes. 



In Figure [2] we show the result for the average current per bond J*(s, L)/(L + 1) as 
a function of s at a = .5, (3 = 2/3 (LD-MC transition line). In this figure, the space-time 
transition at s — 0, rounded by finite size effects, is clearly visible. Qualitatively similar 
behaviour is found for other values of a and (3. 

We next propose a scaling form for the cumulant generating function. That such a 
scaling exists can be expected from analogy with the ring case, equation (20). Moreover, 
the similarity between the CGF and the equilibrium free energy leads us to expect 
that such a scaling can hold near the nonequilibrium phase transitions in the TASEP. 
In writing down a scaling relation for the CGF a natural variable will therefore be 
Aa = a — 1/2, the distance from the transition between the LD and the MC phase (or 
equivalently A/3 = /3 — 1/2 for the HD to MC transition) [35J. These considerations lead 
us to make the finite size scaling ansatz 

//(s, L, Aa) = s(L + I) J* + L~ z H(sL y % AaL Va ) (21) 

where if is a scaling function and y s and y a are two critical exponents that we will 
determine below. The factor L + 1 in the first term on the right hand side equals 
the number of bonds through which the total current passes (including incoming and 
outgoing particles). It replaces the factor L in the ring case. The appearance of the 
dynamical exponent z in the prefactor L~ z is a consequence of the time intensivity of 
fi*. It replaces the space-dimension d in the factor L~ d appearing in the scaling of the 
equilibrium free energy. 
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Figure 2. (Color online) Average current per bond as a function of s at a = 0.50, f3 — 
2/3 for L = 10,20,30,40,50 and 60 (increasing (decreasing) from bottom to top for s 
negative (positive)). 



From (|21| it follows that the average current at s = scales as 

r{s = 0, L, Aa) = {L + l)r + L- z+ys H'{AaL Va ) (22) 

where H' is another scaling function. We can compare this scaling prediction with the 
exact results for large L. Firstly, at the LD-MC transition one has pj5] 



r(s = 0,L,Aa = 0) = ^(1 + ^ + 



(23) 



Comparing with (22), one obtains y s 



z = 3/2. To also determine y a we take the 
derivative of (22) with respect to Aa at a = 1/2. We get 



da 



[s = 0,L,Aa = 0) ~ L Va 



(24) 

one can easily find that 

(25) 



From the exact asymptotic results in the LD-phase 
asymptotically in L 

d £(s = 0,L,Aa = 0) = y^L^(l + ...] 

so that we conclude y a = 1/2. Inserting the exponent values into (22), we obtain the 
final form for the scaling of the CGF near the LD-MC transition. 

(i*(s, L, Aa) = s(L + 1) J* + L- 3 / 2 H(sL 3 / 2 , AaL 1 ' 2 ) (26) 

We have performed several tests of this ansatz with the DMRG-approach and where 
possible using exact results. 
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Figure 3. (Color online) Scaling of J*(s = 0,L, Act) — (L 
AaL 1 / 2 at (3 = 2/3. 



1) J* as a function of 



We start with the scaling for the current at s = as given in (22) since it can be 



checked using the results in [19]. In Figure [3] we present data obtained by numerical 
evaluation of the exact expression for the current for L up to 200. As can be seen the 
scaling is well satisfied, especially from L w 80 on. 

Next we turn to the scaling of the CGF itself. Figure [4] shows a scaling plot of 
(//*(s, L, 0) — s(L + 1)/4)L 3 / 2 as a function of x = sL 3 ^ 2 . These data are determined 
from the DMRG. We used the data for small s and the largest L values that we could 
obtain (L = 60). The curves for various L- values are seen to come closer together with 
increasing L. Unfortunately, the L- values that can be studied with the DMRG at s ^ 
seem to be just outside the scaling region. 

In Figure[5]we present the scaling plot for J*(s, L) as a function of s at Aa = from 
the DMRG data. There is a good collapse of the data, especially at negative s-values. 

Finally, we turn to the variance of the current, which in this non-equilibrium context 
plays a role similar to the specific heat or the susceptibility in equilibrium. How does 



the variance of the current behave near the LD-MC transition? From (26), we derive 
the scaling form at s = 

AV(0, L, Aa) = L 3/2 H 2 (AaL 1/2 ) (27) 

with H2 a scaling function. An exact formula for the variance of the current was derived 
in [20] using the matrix product technique. Due to the complexity of this formula, closed 
expressions for the variance could only be derived at the point a — f3 — 1 in the MC- 
phase, and along the line a + /3 = 1 (LD and HD phase). For these situations it was 
found that A* J grows respectively as L 3 / 2 and as L 2 . We calculated the variance of the 
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Figure 4. (Color online) Scaling of (fj,*(s,L,0) — s(L + 1)/4)L 3 / 2 as a function of 
sL 3 / 2 at a = 1/2,(3 = 2/3 (for L = 50,52,54,56,58 and 60 (increasing (decreasing) 
from top to bottom for positive (negative) s). 
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Figure 5. (Color online) Scaling of J* (s, L, Aa = 0) — (L + 1)J* as a function of 
sL 3 / 2 (for L = 50, 52, 54, 56, 58 and 60 increasing from top to bottom). 



current in several points in the phase diagram using the DMRG technique and verified 
that these results are universal for each phase. Moreover, at the LD-MC transition line 
we found that A*J(0,L, Aa) ~ £i-5o±.02^ cons i s tent with (27) [21]. To be in agreement 
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Figure 6. (Color online) Scaling of the variance of the current at j3 = 2/3. The 
numerical data have been collapsed with the algoritm of [35] . The best collapse shown 
here is found for y a = .57 and z — 2y s = 1.46. 



with the scalings in the LD and MC phase just discussed, H 2 (x) should be constant for 
x 3> and H2(x) ~ x for x 0. In Figure|6]we show our results for A* J(0, L, Aa)L z ~ 2ys 
as a function of AaL Va . Here we used the algorithm presented in [38J to collapse the 
data for the largest L-values. As the figure shows, this is possible but with exponents 
that deviate ten percent from the conjectured ones. This is again because our data are 
not yet in the asymptotic regime where scaling is expected to hold. The prediction of 
linear behaviour for H2(x) when x is negative enough is however clearly verified. 

From (26) one can also obtain a scaling form for the variance of the current as a 
function of s. We have checked that also this form is consistent with the DMRG-results 
as already shown in Figure 3 of [21]. 

In conclusion, the scaling assumption ( 26 ) allows a consistent description of all exact 
and numerical data on the average current and its fluctuations. Due to the symmetries 
of the TASEP we expect a completely similar scaling to hold near the HD-MC transition. 



6. Time dependent current fluctuations 

We first discuss our results for the gap G(s, L) near the LD-MC transition and then 
investigate their implications for the time-dependent current flucutations. 

The behaviour of the gap for s = was first determined numerically using the 
DMRG approach to stochastic operators [33J. Later exact results were obtained from 
an analysis of the Bethe ansatz equations [351 EZ] ■ In the MC phase, it was found that 
G(0, L) goes to zero as L _3//2 , consistent with the criticality of the phase. In the LD and 
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Figure 7. (Color online) Gap G(s, L) as a function of s for different L- values and for 
q = 1/2, /? = 2/3. 



HD phase on the other hand, the gap goes to a constant indicating a finite relaxation 
time. Here we study for the first time the gap as a function of s. 

In Figure [7J we show as a typical result for different L- values at a = 1/2, (3 = 2/3. 
As can be seen the gap approaches \a — /3| = 1/6, the value predicted for large s in 
section 3, already for rather small s- values. 

We have found that all the exact and numerical results for the gap can again be 



described by a scaling function. In analogy with (26) we propose the form 
G(s,L, Aa) = L-^Fistfl^AaL 1 ' 2 ) 



(281 



where F is a scaling function. This form can describe the exact results at s — if -F(0, y) 
goes to a constant for y ^> and goes as y 3 for ?/ C 0. Figure [8] shows our numerical 
results. At s = reliable values for the gap can be obtained up to L fa 100. The scaling 
is already satisfied very well for L fa 50 and the data also show the expected behaviour 
for large \y\. 

In Figure [9j we show the scaling of the gap at Aa = as a function of s. There is 
again a good collapse, especially for negative s-values. 

The scalings that we have found for //* and G imply a scaling for the late time 



behaviour of the time-dependent current fluctuations. Inserting (26) and (28) into (14) 
we obtain 



-sJ(t) 



) = A x e 



i(L+l)J*t 



e L-' i / 2 tH{sL 3 / 2 AaL 1 / 2 ) ^ + ^ e -L- 3 / 2 tF{sL 3 / 2 AaL 1 / 2 ) _^ 



The whole term within square brackets depends on time through the combination tL 3 / 2 , 
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Figure 8. (Color online) Scaled gap G(0, Aa, L)L Z I 2 as a function of AaL 1 / 2 . 
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Figure 9. (Color online) Scaled gap G(s, L)L 3 ^ 2 at the LD-MC transition as a 
function of sL 3 / 2 . 



so that we obtain 

//(s, L, t) = s(L + 1) J* + L~ 3/2 K(sL 3/2 , AaL 1/2 ,tL- 3/2 ) (29) 

where X is once more a scaling function. This equation gives insight in how the time- 
dependent current and its fluctuations reach their asymptotic NESS value. For the 



Current fluctuations in the TASEP 



15 



5 
+ 





1000 



10 



15 



tl_ 



■3/2 



L = 


10 


+ 


L = 


20 




L = 


30 




L = 


40 


□ 


L = 


50 


■ 


L = 


60 


O 


L = 


70 


• 


L = 


80 




L = 


90 




L = 


100 





20 



25 



Figure 10. (Color online) Scaling plot of J(t, L, Aa = 0) — [L + l)/4 as a function 
of tL~^ 2 . The main figure shows the result for intermediate values of the scaling 
variable. The inset shows the data for large t. 



average current at the LD-MC transition, we get for example 
J(t,L, Aa = 0) = (L + l)/4 + Ki(tL~ 3 ^ 2 ) 



(30) 



where K\ is another scaling function. Comparison with (23) gives K\{x) — > 1/8 for 
x — > oo. For the variance of the current one finds 

AJ(t,L,Aa) = L 3/2 K 2 (AaL 1/2 , tL~ 3/2 ) (31) 

The time dependent current and its variance cannot be obtained directly from our 
DMRG-method. Here instead, we have verified relations like (30) and (31) by 
simulations based on the Gillespie algorithm for L up to 100. In Figure (10) we show 
our data for J(t,L,Aa = 0) — (L + l)/4 as a function of tL~ 3 ^ 2 . The initial condition 
was an empty lattice. We simulated 10 4 histories for L < 50. For larger system sizes 
the number of realisations decreased up to 2000 for L = 100. As can be seen the scaling 
is well satisfied. At large times, the current reaches its L-dependent NESS value, which 
for L — )• oo should reach the asymptotic value 1/8 (inset). 

A similar picture arises for the time-dependence of the variance of the current 
though the data are somewhat noisier for the number of histories that we could simulate 
(see Fig. 11). Within their accuracy they are consistent with the prediction (31). 
The variance of the current in the NESS cannot be obtained very precisely from these 
simulations. While consistent with the DMRG-values, the latter are much more precise. 
This is another advantage of the DMRG in comparison with simulations. 
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Figure 11. (Color online) Scaling plot of AJ(t,L, Aa = 0)L 3 / 2 as a function of 
tL-^l 2 at a = 1/2, ft = 2/3. The main figure shows the result for intermediate values 
of the scaling variable. The inset shows the data for t large. 



7. Conclusions 



In this paper we have studied the current fluctuations of the TASEP with open 
boundaries, both in the NESS and as a function of time. 

We have shown that, in the NESS and in the thermodynamic limit, the cumulant 
generating function must be non-analytical at s = implying the existence of a space- 
time phase transition in the TASEP. 

We have proposed the finite size scaling form (29) for the (time-dependent) 
cumulant generating function. Important variables herein are the distance Aa from the 
non-equilibrium phase transition line and s. The form proposed is a natural extension 
of that near equilibrium phase transitions and that found for the TASEP on the ring. 

We have verified our ansatz using available exact results and numerical approaches 
based on the DMRG and the Gillespie algorithm. 

A scaling function for the current generating function has also been derived recently 
for the symmetric simple exclusion process (SSEP) on a ring [39]. The form for the 
CGF derived in that work shows some similarity with our ansatz (21) when, as should 
be expected for a diffusive system, the dynamical exponent z is changed to 2. It would 
be of interest to investigate how general our scaling form is. 

The TASEP was originally introduced as a model for mRNA-translation [12]. We 
have recently shown how fluctuations in the number of proteins produced by one mRNA 
can be related to the time-dependent current fluctuations in a model for translation 
[13]. From the scaling proposed here, it is then possible to derive how, for example, the 



Current fluctuations in the TASEP 



17 



variance in the number of proteins produced depends on mRNA-length. This opens a 
road to a possible experimental verification of our scaling ansatz. 

In this paper we have only studied the two largest eigenvalues of the generalised 
generator H(s). The DMRG also allows the calculation of the associated left- and right 
eigenvectors. From these it is possible to obtain expectation values, such as the density 
at a given site [25J. By tuning the parameter s, it is then possible to see how the typ- 
ical density profile changes for a current that deviates from its average value. We are 
currently calculating these properties. The results will be published elsewhere. 
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